STABILITY OF THE RECOVERY OF MISSING 
SAMPLES IN DERIVATIVE OVERSAMPLING 
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Abstract. This paper deals with the problem of reconstructing a band-limited 
signal when a finite subset of its samples and of its derivative are missing. The 
technique used, due to P.J.S.G. Ferreira, is based on the use of a particular 
frame for band-limited functions and the relative oversampling formulas. We 
study the eigenvalues of the matrices arising in the procedure of recovering the 
lost samples, finding estimates of their eigenvalues and their dependence on 
the oversampling parameter and on the number of missing samples. When the 
missing samples are consecutive, the problem may become very ill-conditioned. 
We present a numerical procedure to overcome this difficulty, also in presence 
of noisy data, by using Tikhonov regularization techniques. 



1. Introduction 

This paper deals with the problem of recovering missing samples of band-limited 
functions and of its derivative via frame reconstruction sampling formulas. A band- 
limited signal is a function which belongs to the space of functions in L 2 (M.) 
whose Fourier transforms have support in [— Functions in this space can be 
expanded in terms of the orthonormal basis of translates of the sine function. The 
coefficients of the expansion are the samples of the function at a uniform grid on R, 
with "density" uj/tt (Nyquist density). Such expansions, called sampling formulas, 
have been generalized by replacing the orthonormal basis with more general fami- 
lies, like Riesz bases and frames, formed by the translates of one or more functions 
(generators). Frames, unlike Riesz basis, are overcomplete and their redundancy 
provides a perfect tool for the recovery of missing data. 

The problem of recovering missing samples has been investigated first by P.J.S.G. 
Ferreira in (FJ, where it is shown that, in the case of a generalized Kramer sampling, 
under suitable oversampling assumptions, any finite set of missing samples can be 
recovered from the others by solving a linear system (I — S)x — b, where the 
matrix S is positive definite. Moreover the author studies the eigenvalues of the 
matrix S in dependence of the oversampling parameter and the number of missing 
samples. Successively D.M.S. Santos and Ferreira considered the case of a two- 
channel derivative oversampling formula obtained by projecting the generators of 
a Riesz basis of the space B u and their duals into the space B bJa with cu a < lo [SF\ . 
With this technique they obtain a pair of dual frames, although the dual frame 
is not the canonical one. In their paper the authors show that a finite number of 
missing samples either of the function or of its derivative can be recovered, solving 
in each case a non singular linear system similar to the one-channel case. However 
the authors do not consider the case when samples both of the function and the 
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derivative are simultaneously missing. 

The technique of Santos and Ferreira has been applied to more general two-channels 
by J. M. Kim and K. H. Kwon |KK| , who gave sufficient conditions for the recovery 
of missing samples, from a single channel and from both channels. However, as the 
authors observe, these conditions do not apply to the derivative sampling formulas 
studied in jSFj . 

Successively, in |DPj . the author gave new derivative oversampling formulas of 
any order, where the dual generator is the canonical dual and in [DPI] found an ex- 
pression of the dual in the two-channel case. The use of the canonical dual allowed 
the author to prove in 1DP1] that simultaneously missing samples of the function 
and of the derivative can be recovered. However the generators of the canonical 
duals are much less explicit than the non-canonical dual used by Santos and Fer- 
reira see [DTT1 (5.2) p. 178] and plO| . 

In this paper we study the problem of the recovery of missing samples of the function 
and of its derivative in Ferreira's two-channel derivative formula. The technique 
we use, first proposed in |SFj . extends in a natural way that for one channel in [F] 
and consists in solving a system (/ — S)X = C, where the matrix S (see (3.8)) is a 
block matrix depending on the dual generators and on the position of the missing 
samples, and the unknowns X are the missing samples of the function and its deriv- 
ative. To the best of our knowledge, no results are known about the possibility of 
solving this system, i.e. of recovering simultaneous missing samples of the function 
and of its derivative. 

We obtain estimates of the minimum and maximum eigenvalues of the block sub- 
matrices of the matrix S, show that in some cases the matrix reduces to a lower 
triangular matrix and provide several numerical experiments on the dependence of 
the eigenvalues on the parameters of the problem, like the oversampling parame- 
ter r and the number of missing samples. Moreover, we experiment the recovery 
of missing samples and the reconstruction of a signal, paying particular attention 
to the computational aspects, analyzing the ill-conditioned problem of contiguous 
samples, also with noisy data. 

The paper is organized as follows. In Section 2 we establish notation and collect 
some results to be used later. In Section 3 we describe the technique for the recovery 
of the missing samples. Section 4 is dedicated to the study of the eigenvalues of the 
coefficient matrices in the system, mostly when the missing samples are equidistant. 
In this case, we find estimates of the minimum and maximum eigenvalues of the 
diagonal submatrices Su and S22 of S (see (3.8)). Moreover, we find that in some 
cases the matrix S is lower triangular and that half of its eigenvalues are equal to 
2r — r 2 , half to r 2 , where r is the oversampling parameter. This extends to two 
channels a result in [F]. We also present some numerical experiments supporting 
the theory. Finally, in Section 5, we analyze the case of contiguous missing samples 
both for one and two channels, when the problem may become very ill-conditioned. 
We present a numerical procedure to solve it, also in presence of noisy data, via 
regularization techniques. 
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2. Preliminaries 

In this section we collect some results on frames to be used later. We begin by 
introducing some notation. The Fourier transform of a function / in L 1 (R) is 

^7(0 = /(0 = ^ / f(t)e-^dt. 

In this paper vectors in C N are to be considered as column- vectors, however, to 
save space, we shall write x = (xi 1 x 2 , ■ ■ ■ ,x N ) to denote the column- vector whose 
components are Xi, . . . , Xn. Let t Q be a positive number; we shall say that a 
subspace H of L 2 (M.) is t Q -shift-invariant if H is invariant under translation by t a . 
Given a subset $ = {ipj, j = 1, 2} of H, we denote by E$. to the set 

E<s>,t a = {TntoVj neZ, j = 1, 2}, 
where r a f{x) = f(x + a). The family -E$,t is a frame for H if and only if there 
exist two constants < A < B such that A||/|| 2 < Y^ij=i Snez K/> ^tafj)] 2 < 
P||/|| 2 for all f £ H; the constants A and B are called frame bounds. Denote by 
T|_ to : if -+ ^(Z;^) the adjoint of T* >to , defined by (T| to /),(n) = (f,T nto <pj), 
n e Z, j = 1, 2. The operator T $ t T| s : PT — > PT is called frame operator. Denote 
by $* the family {<p*,j = 1,2}, where 

(2-1) = (T* ito Ti it J-Vi .7 = 1,2. 

If -E*,t is a frame for then P$*,t D is also a frame, called the canonical dual 
frame, and T^ ^T^, ^ = T^, t Tl to = I; explicitly 

2 2 

(2.2) / = x;j;{/,r nto ^T nt ^ j ^Ev-^ft^^i v ^ e 

.7=1 n£Z j = l n£Z 

The elements of $ are called generators of the frame and the elements of canon- 
ical dual generators. Given a frame P$,t of a Hilbert space H, the expansion of an 
element in terms of the generators is not unique; a frame {Tkt Vii T kt ( P2 *= ^} 
such that 

2 

(2-3) / = r ntoV / )r nto ^ V/ e H 

]=1 neZ 

is called dual frame for E$ tto . We remind that the coefficients {(/, T nto tp*),n,j G Z} 
with respect to the canonical dual frame have the minimal ^ 2 norm among all the 
sequences that represent an element / in terms of a given frame P$,t - If the family 
E<p ito is a frame for H and the operator T $ t is injective, then E^, to is called a 
Riesz basis. 

Let P$,t be a frame of a i -shift-invariant Hilbert space -ff with canonical dual 
generators $*. Denote by P the orthogonal projection on a t Q -shift-invariant sub- 
space V of H, then Pp$,t and E P ^* jto are families of dual frames for V, i.e. 

2 2 

/ = E Z</> Tnt PV*)T nto Pipj = £ T nto P<P j )T nto P<p t j V/ € V. 

j=l neZ j = l neZ 

Note that -Ep$*,t is not necessarily the canonical dual. 

In this paper we study families P$,t , $ = {<£i, </?2} where t/?i and if 2 have Fourier 
transform <fi(x) — X[-u,u]i x ) an d ^2( x ) = i x X[~u>.ui] ( x ) an d *o an d w are such that 
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u) < y~ < 2uj. The interest of these families in applications resides in their con- 
nection with derivative sampling formulas for the space B^. To simplify notation, 
throughout the paper we shall set 

to 



The family Ey >to , <3> = {^1,^2} defined by 
(2.4) Mx) = X[-h,h]( x ) 



is a Riesz basis for the space Bh (see 
dual generators are 

(2.5) rx{x) = \{l-~)X[-h,h]( x ) 
hence the dual generators are 



ip 2 (x) = ixx[-h,h]( x ) 
p. 135]). The Fourier transforms of the 



1 . 2 - hx 



sine ( — ) 
2tt V 2 1 



ip 2 (x) 



1 2/h x \ 
-=xsmc ( — ), 

/Z7T ^ 



where sinc(x) = sin(x)/x. For any f € Bh the coefficients of the expansion in (2.2) 
are 



(/, T- nfo ipi) = \f2nfinto) 



(f, T-ntM 



thus the expansion formula (2.2) becomes 



(2.6) f(x) = (/(nt )sinc 2 ^(x- nt a ) + f'(kt )(x - nt a ) sine 2 ^ (x - nt ) ). 

n 

for any / G .B^. Note that the convergence is uniform. This formula is called a 
first order derivative sampling formula. The sampling frequency for an A^-channcl 
formula is lo/Ntt and is called Nyquist frequency. By using a frame instead of a 
Riesz basis, in |SF] the authors obtained a two-channel derivative formula where 
the sampling frequency in each channel is greater than the Nyquist frequency to /2tt 
(first order derivative oversampling formula). The frame is obtained by projecting 
on J5 U the Riesz basis E^ to for B^ where h — 2ir/t and to < h. We denote by r 
the ratio 



(2.7) 



to 
h 



t0t o 

27 



Since the projection commutes with the translations, to project the family E^, 
and the dual family, it is sufficient to project on B u the generators and the dual 
generators. Thus, by (2.4), we obtain that the generators of the frame E$ t t are 
defined by 



(2.8) 



ip 2 (x) = ixxi- u , u ](x) 



and by (2.5) we obtain the Fourier transform of the dual generators 
_ _ ^2 



(2.9) &(x) = -(1 - Hx|) X[ - 

tO LO 



,](x)(x) tp 2 (x) = z-^sign(x)x[_ w ,w](x) 
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A simple calculation shows that 

1 



(2.10) 



ipi(x) 



. \ . / \ 2 • 2 \ 

2r (1 — r)sinc(w3;) + r sine ( — ) 



I 2 . „ M 

xr smc — . 
v 2 y 



2tt 



Thus we obtain the following derivative oversampling formula 

(2.11) f(x) = V2^J2 (f( kt o)Tkt $i ~ f'(kt )T kto $ 2 ) V/ e B w 

[SFj . We observe that 1/r = 2ir/ujt is the ratio between the sampling frequency 
and the Nyquist frequency and that r£ (0, 1). We shall be mainly interested to the 
case 1/2 < r < 1, since for 0<r<l/2itis possible to use one channel separately 
for the function and for its derivative }SF] , Note that if r is close to 1, the frame is 
close to a Riesz basis, while if r is small the frame is very redundant. 



3. The system for the recovery of missing samples 

In this section we briefly describe the method for the recovery of a finite num- 
ber of missing samples via two-channel derivative oversampling. We may rewrite 



equation (2.11 ) as 

2 

(3.1) f{x) = V27rJ2 ^(-ly^f^HktoWiix - kt ) 

i=l kez 

and computing the derivative of both sides, we obtain 

2 

(3.2) f(x) = V2^J2 ^(-l)*- 1 /'"'^)^ " kt ). 



i=l k£Z 



Here 
(3.3) 



2r 1 



2irx 



(1 — r)( cos(o;cc) — sincere)) + 

+ rsinc(-x)(cos( — ) - smc( — )) 



(3.4) &{x) 



1 o . ujx r 
— r sine — ) 2cos( — 
v 2 L v 2 ' 



2tt 



smc( — )) 



Observe that <p' 2 (0) = 0. 

Let U = {h,k.-.,lN} C Z and let {f(ljt ), f'(ljt ) I < .j^ N} be the corre- 



sponding set of missing samples. By evaluating equations (3.1) and (3.2 1 in l]~t 



and separating the unknown samples from the known ones, we obtain 

2 N 

f(i k t Q ) - v^^^rt-ir- 1 /^ 1 ^) (piiih-i^to) = Ck 

(3.5) 1=13=1 

v ' IN 

fUkto) - v^^r^-ir- 1 /^ 1 ^)^'^ - Wo) = cfc+jv, 



i=i j=i 
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1 < k < N, where 



2 



(3.6) 



i=l n& 
2 



Ck+N 



i=l n£i 



Equations (3.5) form a system of 2N equations in the 2N unknowns 

(3.7) f(ht ), . . . , f(l N t ), f{ht ), . . . , f{l N t ) 

which can be written in matrix form: denote by S — S(U, t) the real matrix 

S21 S22 

where Su, Si 2, S21, S22 are the submatrices whose entries are 

Sn(k,j) = V2tt tpi((lk - lj)t a ) S 12 {k,j) = -V2ir (f 2 ((h - lj)to) 

(3.9) 



(3.8) 



S 



S2i{k,j) = V2n tpi'((lk - lj)t ) S 2 2(k,j) 



27T ip 2 ((Ik ~ lj)to) 



k,j = 1, 
(3.10) 



. , N. Then the system ( |3.5| may be written 
(I-S)Z = C 



where / is the 27V x 2N identity matrix, C = (cfc)^j is defined by (3.6) and 
Z = [Xi, X 2 , ■ ■ ■ Xn, Yi, Y 2 , . . . Yat) . The unknowns Xj, are the missing samples of 
the function and the are the missing samples of the derivative. 
Notice that, by (2.10), (3.3) and (|3.4|, if r = 1 then S is the identity matrix. 



The four submatrices are real, Su and 5*22 are symmetric, while S12 and S21 are 
antisymmetric. Moreover, if the distance between two consecutive Ik is constant, 
then the four submatrices are Toeplitz, but S is not Toeplitz. 
Suppose that in (3.7) only the values {Xi = f(lit ), i = 1, . . . , N} are missing, 
while all the values of /' are known. By considering only the first N equations 



of the system (3.10) and by separating the known from the unknown samples, we 
obtain 

(7 - S U )X = C x + S 12 Y, 



where / is the N x N identity matrix, X = Y = (/'(/&£„)) and 

Ci = {ck)k=i- This equation may be rewritten 

(3.11) (I-S U )X = B 1 

where Bx = C\ + S\2 Y . Similarly, if in (3.7) only the values {Yj = f'(lit ), i = 
1, . . . , N} are missing, and the samples of / are known, one solves the system 

(3.12) (I~S 2 2)Y = B 2 

in the unknowns Y i: i — 1, . . . , N , where C2 = {ckjk=i + N and B 2 = C2 + S21X. 
If r g (0, 1), then the eigenvalues of the symmetric matrices Su and ^22 are in the 
interval (0, 1) SF . This implies that the recovery of a finite number of missing 
samples of the function is possible, if all the samples of the derivative are known 
and, vice- versa, samples of the derivative can be recovered, if all the samples of the 
function are known. In the next section we shall find bounds for the minimum and 
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maximum eigenvalues of the matrices Su and S^- 

To the best of our knowledge no results are known about the possibility of solving 
system (3.10), i.e. of recovering simultaneous missing samples of the function and 
of its derivative. Notice that the system is solvable if and only if 1 is not an 
eigenvalue of the matrix S. On the basis of numerical evidence we conjecture that 
all its eigenvalues are real and that they lie in the interval (0, 1) for all r £ (0, 1). 
In Section 4 we present several numerical experiments supporting our conjecture. 



4. The stability of the matrices. 
This section is dedicated to the study of the eigenvalues of the matrices S, Sn and 



S22 in the two-channel system (3.10). First, we briefly summarize the one-channel 
case and recall some stability results matrix obtained in [FJ. The one-channel 
formula is 

(4.1) f (x)= ^f {nto) Vf£B„ 

where t a < r-. In this case the oversampling parameter is r = ujt /TT. Let U = 
{Zij Zi, . . . , £jv} be a set of integer numbers and 

(4.2) {f(l k t ),k=l,...,N} 

the set of missing samples; then the system to be solved is 

(4.3) (I - R)X = B 
where I is the identity N x N matrix 

(4.4) R(j,k) = rsmc(irr{l j -l k )), j, k = 1, . . . , N, 
X = {f(l k t ),k=l,... > N},B = (b u b2,---,b N ), and 

(4.5) bj = r f(kt )smc(Trr{l 3 - l k )) j = 1, . . . ,N. 

The matrix R is symmetric and positive definite. In [FJ the author shows that all 
its eigenvalues are in (0,1) and observes that there are two situations that lead 
to a ill-conditioned problem: when r is close to 1 and when the integers l k are 
contiguous. In the first case, the frame is close to a Riesz basis, when the recovery 
is impossible and R becomes the identity matrix. In the second case the maximum 
eigenvalue of the matrix R grows rapidly with N; moreover it can be close to 1 also 
for small values of N (see Figure 5 in [FJ). In both cases the spectral condition 
number of the matrix I — R, which controls the propagation of errors on the data, 
becomes very large. 

Next, we investigate the eigenvalues of the matrix S in the two channel case (see 



(3.8)) and of its submatrices Sn and S22. Since the latter are real symmetric 
matrices, their eigenvalues give their condition number. We shall denote by \_a\ the 
maximum integer less than or equal to a. Given a. N x N matrix, we shall denote 
with Xj (A) , j = 1, . . . , N its eigenvalues and by \ m in(A) and \ m ax(A) its minimum 
and maximum eigenvalue. 

Let U — {hih ■ ■ ■ Jn} C Z and let {f(ljt ), f'(ljt ) 1 < j < N} a set of missing 
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samples. Then 



X m in(Sll) < 2r — r 2 < X max (Sii) 

X m in{S 22 ) < f 2 < X max {S22) ■ 



where r is the oversampli ng pa rame ter ( 2.7). Indeed, by (3.9) and (2.10), the trace 
of S\i is N(2r — r 2 ); by (3.9) and (3.4) the trace of S 2 2 is Nr 2 . Moreover, since 



the entries of the matrix S are real, its eigenvalues are complex conjugates, thus its 
trace is N(2r 

(4.6) 



) + Nr 2 = J2t=i Re(Xi(S)); hence 
Re(X mm (S)) < r < Re(X 

max 



We observe that, if the oversampling parameter r tends to 1, i.e. the frame tends 
to a Riesz basis, then S tends to the identity matrix. 

Following Ferreira, we shall now consider the case in which the set that locates 
the positions of the missing samples is U = {m i±, raij, . . . , m ijv}, where m is 
an integer and X — i 2l ■ ■ ■ *tv} is a set of integers; in the following we shall 
denote such sets by ml. The interest for studying these cases lies in the technique 
of interleaving the samples of a signal, prior to their transmission or archival; the 
advantage of this procedure is that the transmitted (or stored) information becomes 
less sensitive to the burst errors that typically affect contiguous set of samples ([F]). 
We shall investigate how the stability of the method depends on the interleaving 
factor m. In Proposition |4.1| and Theorem |4.2| below we find estimates for the 
eigenvalues of the matrices Sn,S 22 , thus generalizing a result of Ferreira for one 
channel [0 Theorem 1]. 

First we consider the case mr integer. The following proposition, which extends to 
two channels a result in [F] , shows that in this case the matrix S is lower triangular 



and that N eigenvalues are equal to 2r 



and N are equal to r . Hence if 



is rational r = p/q, one can take m equal to a multiple of q and obtain a lower- 
triangular matrix S. 

Proposition 4.1. Let m be a positive integer; suppose U = {mix, mi^,- ■ • , mi^}, 
where = 1,...,N are integers and let r be a real number in (0, 1). If mr is an 
integer, then 

S u = (2r - r 2 )I S 22 = r 2 I S l2 = 0. 
Moreover the entries of the matrix S 2 \ are 

\ {I - r)u / {mir(ik - ij)) if k^j 
[0 if k = j. 



Proof. Since t a = 2 r Kr/oj, by (3.9) with = mij~, one obtains that 

r 

,— i 



Sn(k.j) = V2TTipi(2nm—(i k ~ i,-)); 



hence by (2.10) since mr is an integer 

S u (k,j) - (2r (1 - r) + r 2 )S jtk l<k<N, l<j<N 
where 6j t h is the Kronecker symbol. This shows that Su = (2r — r 2 )/jy Similarly, 



by using (3.9) and (3.4), we obtain that S; 



r I and Si 2 = 0. From (3.3) one 



finds that the off-diagonal entries of the matrix S21 are 

S 2 i(M = (l-r) k^j. 

mir(i k - ij) 
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The diagonal entries in S21 are equal to zero, since <p'i(0) = 0. □ 

In the following theorem we find estimates for the minimum and maximum 
eigenvalues of the matrix S n and S 2 2- The estimates do not depend on the number 
N of missing samples and are simple to compute. 

Theorem 4.2. Suppose that < r < 1 andU = ml, where I — { i\, 12, ■ ■ • , In} 

and m is a positive integer. If mr is not integer, then for j = 1, . . . ,N 

(4.7) an < X^Su) < fa 

(4.8) a 22 < Aj(5 22 ) < fa 
where 

D D 1 1 

an = — {!--. 5—), fa = au-\ — 

m 4m 4m m 



and D = \ 2mr\ . 



Proof. From (2.9), by Fourier inversion and a change of variable, we obtain 
(4.9) = -±= f (1 - \y\)e^ r dy. 

V ZTT J-r 



Hence by (3.9) with l k = mif., 

Su(k,j)= I (l-\y\)e 2 ™ m ^-^dy. 



Let x — (x , X\, . . . , xn-i) be a column-vector in C N and denote by x* its conjugate 
transpose. Then 

x*S llX = J2 x k Sn(k,j) Xj = / (1-M) XkXje 2mm y^-^dy 

k,j=0 r k,j=0 

(l-\y\)\P(y)\ 2 dy, 

-r 

where P{y) is the 1/m-periodic function 

JV-l 

(4.10) P(y) = x k e- 2 ^ myik . 

fc=0 

By splitting the integral in two parts and changing variable the left hand side of 
the above equation may be written 



(4.11) x*S n x= / (l-y)M(y) dy 

Jo 

where we have set 

(4.12) M(y) = \P(y)\ 2 + \P(-y)\ 2 . 



Next we write the integral in (4.11 ) as a sum of D + 1 integrals 

r(k+l)/2m r-r 

(4.13) x*S n x =Y, (1 - y)M(y) dy + (1 - y)M{y) dy. 

k=0 "' fc / 2m J D/(2m) 
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Since M > we majorize 1 — y in each integral and obtain the estimate 

D-l , f (fe+l)/2m 



(4.14) x*S llX <^{l 



fc=0 



2 m 



M{y) dy + (l 



k/2m 



2m' 



M(y) dy. 



D/(2r. 



Next we prove that for all k = 0, . . . , D — 1 

»(fe+l)/2m 

(4.15) 



1 



M{y) dy = 

k/ '2m 171 



Indeed by ( 4.12[ ), by changing variables and using the (l/?7i)-periodicity of P, 



(fc+l)/2m 



M(y) dy 



k/2r. 



(fc+l)/2m 

k/2m 
(fe+l)/2m 

k/2m 
(fe+2)/2m 

k/2m 



\P(y)\ 2 dy 



\P(y)\ 2 dy 



\P(y)\ 2 dy 



-k/2r, 

-(k+l)/2m 
(fe+2)/2m 

(fc+l)/2m 
1 / m 



\P{V)\ 2 dy 



\P(V)\ 2 dy 



\P{y)V dy. 



By (4.10) and the Planchcrel formula, since x has norm equal to 1, we have 

r-l/m 



(4.16) 



\P(y)V dy 



1 



This concludes the proof of equation ( 4.15[ ) . Next we prove that 

1 



(4.17) 



M{y) dy < 

D/(2m) m 



Indeed, from (4.12), by changing variable in the integral and using the 1/m- 
periodicity of P, we obtain 



D/(2r. 



M(y) dy = 



D/(2r, 



D/(2r, 



\P(y)\ 2 dy 



\P(y)\ 2 dy 



D/(2r. 
D/(2m) 

D/r 



\P(-V)\ 2 dy 



\P(y)\ 2 dy 



\P(y)\ 2 dy. 

' D/m—r 

The length 2r — D/m of the interval of integration is less than the period 1 jm of 
P, thus inequality (4.17) follows from (4.16). 
From ( |4.14[ ), by using ( |4.15| ) and ( |4.17 ), we obtain 

D-l 



(4-18) x*S 11 x<-\J2(1- 7 ^) +(!-£)] =5(1- D 
m L ^— ' 2m 



k=0 



2m 



1 ) + 1 

4m Am m 



Thus we have proved the upper bound for the eigenvalues of S\\. From (4.13), by 
observing that the second integral is positive and using (4.15) we obtain 



D-l „(fe + l)/2m -, D-l 

x*S n x > / (l-v)M(y) dy > - V (l 

k=o Jk / 2m m k=0 



fc + F _ D, D___L_. 
2m m 4m 4m 



This proves the lower bound and concludes the proof of inequality (4.7). The omit 
the proof of inequality (4.8) which is similar. □ 
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r 


an 


Amm(S'll) 


Knax(Sll) 


Pn 


0.55 


0.719 


0.768 


0.811 


0.844 


.6 


0.773 


0.813 


0.859 


0.898 


.7 


0.859 


0.903 


0.926 


0.984 


.8 


0.891 


0.946 


0.967 


1.016 


.9 


0.929 


0.984 


0.998 


1.055 


.95 


0.936 


0.996 


0.999 


1.063 



Table 1. A min (5 , n), \ ma x(Sn) and their estimates for U = 
{0, 8, 16, 24}, m — 8 and several values of r. 



r 


«22 


Amm (S22) 


Amax(*5'22) 


/?22 


.55 


0.219 


0.271 


0.315 


0.350 


.6 


0.281 


0.317 


0.391 


0.427 


.7 


0.430 


0.470 


0.535 


0.603 


.8 


0.516 


0.594 


0.659 


0.709 


.9 


0.711 


0.766 


0.871 


0.932 


.95 


0.820 


0.877 


0.962 


1.056 



Table 2. A min (5 , 22), \max{S22) and their estimates for U = 
{0, 8, 16, 24} and several values of r. 



In Table [T] and Table [2] we compare the minimum and maximum eigenvalues 
of the matrices Su and S22 with their estimates given by Theorem |4.2| for U = 



{0, 8, 16, 24} m = 8 and various values of r. From Proposition 4.1 one can see that, 
if mr is integer, then the eigenvalues of £22 are smaller than the eigenvalues of 
S11. The following corollary shows that this is also true when mr is not an integer, 
provided that m is sufficiently large. 

Corollary 4.3. Let m be a positive integer, U = ml and X = ii-, • • • , «Af}, 
< r < 1. If mr is not an integer and m > (1 + 2r)/(2r(l — r)), then 

(4.19) A 

lilii.r (S22) < A m j n (S'ii). 



Proof. By gj]) and (gTSj) it is sufficient to prove that if m > (1 + 2r)/(2r(l - r)), 
then P22 < ocn . This inequality is equivalent to D 2 /2m + D/2m + r — D < 0. Since 
2mr — 1 < D < 2mr, we have 

D 2 D f 

V Yr - D < 2mr(r - 1) + 2r + 1 

2m 2m 

from which the corollary follows. □ 

Numerical experiments show that for m < (1 + 2r)/(2r(l — r) the maximum 
eigenvalue of S22 may be larger that the minimum eigenvalue of S11. 

We shall now describe the behavior of the eigenvalues of the matrix S in de- 
pendence of the parameters r and N. The numerous experiments that we have 
performed suggest the conjecture that the eigenvalues of this matrix are real, pos- 
itive, and less than 1. In what follows we have ordered the eigenvalues so that 
A, < Aj + i for i = 0, . . . , N — 2. Figure [T] shows the largest eigenvalue of S as func- 
tion of the number N of contiguous points, for various values of r. Note that, as r 
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Figure 1. Largest eigenvalues of S as function of TV for X 
{0, 1, ■ ■ • ,N — 1}, for various values of r. 



r 


cond(J - S) 


r 


cond(7 - S) 


.1 


8.571 c+01 


.6 


2.829 e+13 


.2 


5.870 e+04 


.7 


1.661 e+16 


.3 


6.187 e+05 


.8 


2.474 e+17 


.4 


1.133 e+08 


.9 


4.096 e+17 


.5 


3.513e+10 


1 





Table 3. Condition number of I 
0.1,0.2,...,!. 



S for U = {0,1,2, ...,9} and r 



approaches 1, the largest eigenvalue gets close to 1, i.e. the smallest eigenvalue of 
the matrix I — S tends to zero. This happens even for a small number of missing 
samples, as for the one-channel case (see [F], Figure 5). Thus, since the spectral 
condition number of a matrix is greater or equal to the ratio of the absolute values 
of the largest and the smallest eigenvalues, cond(I — S) becomes very large. In 
Table [3] we show the spectral condition numbers of / — S for U = {0, 1, 2, • • • , 9}, 
for various values of r. 

Figure [2] shows some of the eigenvalues of S, namely Xi,i — 1,5,10,11,15,20, as 
functions of the parameter r for 2 = {0, 1, 2, 3, • • • , 9}, m = 6 (left) and m = 10 
(right). The dashed lines are the parabolas 2r — r 2 and r 2 . One can see that for 
larger values of m the graphs of the eigenvalues concentrate around the graphs of 
the two parabolas. We observe that in the one-channel case the eigenvalues have a 
similar behavior (see Figure 7 in [F]). 

We shall now discuss the behavior of the eigenvalues of the matrix S as the pa- 
rameter m grows. By using (3.9) with I). = mi)., k = 1,.,.,N, (2.10), (3.3) and 



(3.4), one can see that, when the parameter m tends to infinity, the off-diagonal 
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0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 



0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 09 1 



Figure 2. Six eigenvalues of S as functions of r for U = 
{0, 1,2,3, . . . ,9}, m = 4 (left) and m = 10 (right). Dashed lines: 
parabolas r 2 and 2r — r 2 




Figure 3. Some of the eigenvalues of S as a function of m for X 
{0,1,. ..,9}. r = .7 



entries tend to zero; for each e > denote by m e the positive integer (depending 
on N) such that, if m > m e , the off-diagonal entries of the matrix are less than e. 
For such m e , by the Gershgorin theorem, N eigenvalues lie in the circle of center 
2r — r 2 and radius e and N eigenvalues lie in the circle of center r 2 and radius e. 
Figure [3] illustrates this behavior: it shows the eigenvalues A, i — 1, 5, 10, 11, 15, 20 
in dependence of the parameter m for r — .7 and X — {0, 1, . . . , 9}. Note that the 
limits of the eigenvalues are 2r — r 2 and r 2 . 



5. Reconstruction: numerical results 

In this section we present some numerical experiments on the recovery of miss- 
ing samples and on the reconstruction of a signal via oversampling formulas. In 
particular, we analyze the case of contiguous missing samples, when the problem 
may become very ill-conditioned, depending on the parameter r and the number of 
missing samples. Indeed if r is close to 1, or the number of missing samples is large, 
the condition number of the matrix J — S is large, so that the system amplifies 
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the errors in the data. We solve this problem, also in presence of noisy data, by 
applying a regularization technique typical of the treatment of inverse problems 
(see |BB) . Section 5). These techniques consist in considering a family of approxi- 
mate solutions X\ depending on a non- negative parameter A, called regularization 
parameter. When the data are noise-free, the solution X\ converges to the exact 
solution as the regularization parameter tends to zero. In the case of noisy data, 
one can obtain an optimal approximation of the exact solution for a positive value 
of the parameter. 

Following Ferreira, in all our experiments we shall use the test function 
(5.1) g{x) = (sinc(7r(ai - 2.1))) - 0.7(sinc(7r(a: + 1.7))). 

which has band [— tt,tt] |Fj. 

First we present some numerical experiments for the recovery of contiguous missing 
samples via the one-channel formula (4.1 1. To compute the sum in (4.5), we must 
truncate it to the values k such that \k\ < M, for some integer M, thus introducing 
an error. In [F] the author makes the choice M = 40 to reconstruct the signal 
(5.1) when U = {0, 1, ... ,5} and r = 0.6 (hence t a — 0.6). One can verify that 
for M = 40 the norm of the truncation error is of order 10 -3 and that choosing 
M = 500 this error is reduced to the order of 10~ 4 , while, due to round-off errors, a 
larger value of M does not reduce it further. In Table|4]we show the values Xk = hto 
(column 1) the exact samples g{xk) (column 2), the computed samples g* k and the 
relative errors e* k with M = 40 (columns 3 and 4), the computed samples gk and the 
relative errors with M — 500 (columns 5 and 6). The improved results obtained 
in the latter case can be explained by observing that the condition number of the 
matrix I — R. which controls the propagation of the errors, is 3.07 • 10 4 . We have 



Xk 


g(xk) 


9 k 




9k 


ek 


0.0 


0.1529 


0.1132 


0.2598 


0.1498 


0.0199 


0.6 


-0.2906 


-0.5344 


0.8390 


-0.3096 


0.0653 


1.2 


0.0856 


-0.4833 


6.6498 


0.0410 


0.5206 


1.8 


0.9221 


0.2132 


0.7687 


0.8664 


0.0603 


2.4 


0.8416 


0.3498 


0.5844 


0.8029 


0.0459 


3.0 


0.0710 


-0.0872 


2.2288 


0.0585 


0.1751 



Table 4. Recovered samples with the one-channel formula, r=0.6: ex- 
act samples (col 2), computed samples and relative errors with M = 40 
(col.s 3 and 4) and M = 500 (col.s 5 and 6). 



also plotted the graphs: in Figure [3] we compare the original and the reconstructed 
signal in both cases. For M = 500 we have zoomed the graphs, representing them in 
the interval [—5, 5], to make them distinguishable. By the preceding considerations, 
in all the experiments described below, we have chosen M — 500 to reduce the error 
due to the truncation error. 

Next we consider the problem of data affected by noise; we choose again r = 0.6 
and U = {—2,-1,0,1,2,3}. The condition number of the matrix I — R is still 
3.0769 • 10 4 , sinc e it does not depend on the position of the contiguous missing 
samples (see (4.4)). We have introduced a random noise of order 10 -2 on the data 
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-25 -20 -15 -10 -5 5 10 15 20 25 -5 -4 -3 -2 -1 1 2 3 4 5 



Figure 4. Original (solid line) and reconstructed (dotted line) signal 
g, r = .6, U = {0, 1, ... ,5}. Truncation M = 40 (left) and M = 500 
(right), on the interval [—5, 5]. 



and have solved the system (3.10) by using the Tikhonov regularization technique. 
Since an estimate of the norm of the noise is known, it is possible to use the 
discrepancy principle to choose the value of the regularization parameter A (see 
[BB] and [M])- Table [5] summarizes the results of this experiment; columns 1 and 
2 contain the values Xk = Ikto and the exact values g(xk) of the missing samples, 





g{xk) 


9k 


efe 


9l 


e fc 


-1.2 


-0.5237 


50.9126 


98.2227 


-0.4852 


0.0734 


-0.6 


0.1580 


180.6426 


1142.5116 


0.2141 


0.3550 


0.0 


0.1529 


306.9469 


2006.8097 


0.1430 


0.0652 


0.6 


-0.2906 


307.9712 


1060.8440 


-0.3608 


0.2416 


1.2 


0.0856 


183.2893 


2141.4934 


0.0435 


0.4918 


1.8 


0.9221 


53.8130 


57.3621 


0.9248 


0.0030 



Table 5. Recovered samples with the one-channel formula, r=0.6, with 
noise added: exact samples (col 2), computed samples and relative errors 
without regularization (col.s 3 and 4), computed samples with regular- 
ization and the relative errors (col.s 5 and 6). 



columns 3 and 4 contain the recovered samples gk and the relative errors obtained 
without regularization, while columns 5 and 6 show the recovered samples g 7 k and 
the relative errors ej* obtained via regularization. As in the previous example by 
comparing columns 2 and 5, one can see that the absolute errors are very small: 
the plots of the real and reconstructed signal would be indistinguishable. 

We shall now present some experiments that involve the two-channel formula. 
To compute the right-hand side of (3.10), we have truncated the sums in (3.6) 
retaining only the terms for which \n\ < 500. As in the one-channel case, this leads 
to an error of order 10~ 4 . Our experiments show that, as in the one-channel case, 
the recovery and the reconstruction of missing samples of the function and of its 
derivative is rather efficient if the distance between the missing samples is large 
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and the oversampling parameter r is not too close to 1. To give an example, for 
X = {—1, 0, 1, 2, 3, 4}, m = 4 and r — 0.7, the maximum absolute error is less then 
8 • 10~ 4 . The error increases when the parameter m gets smaller, or r larger. 

Next, we experiment the recovery of consecutive missing samples of g and g' 
when U = {—2,-1,0,1,2,3}, r = 0.6 (hence t a = 1.2). The results are shown 
in Tables [6] and [7| In this case the condition number of the matrix I — S in 
(3.101 is 3.67 • 10 7 , much larger than in the one-channel case. This has the effect 
of propagating strongly the truncation error, thus producing huge errors of the 
solution, much worse than in the one-channel case. Thus, with this distribution of 
missing points and this value of the parameter r the system is very ill-conditioned; 
an application of the Tikhonov regularization technique with discrepancy principle 
leads to the results in columns 5 and 6 in Table [6] and [7J In Figure [5] we show the 
plots of the original and the computed signal g. 





g(xk) 


9k 


e/c 


9l 




-2.4 


-0.1868 


-0.0580 


-0.6893 


-0.0785 


0.5799 


-1.2 


-0.5237 


2.1433 


-5.0929 


-0.4569 


0.1274 


0.0 


0.1528 


10.3983 


67.0174 


0.1832 


0.1989 


1.2 


0.0856 


12.0422 


139.7624 


-0.1897 


3.2170 


2.4 


0.8416 


5.1810 


5.1561 


-0.1480 


1.1758 


3.6 


-0.1782 


0.1425 


-1.7997 


-0.4405 


1.4720 



Table 6. Recovered samples of g with two channels, r = 0.6, no noise 
added: exact samples (col. 2), computed samples and relative errors 
(col.s 3 and 4), computed samples and the relative errors (col.s 5 and 6) 
with regularization. 



Xk 


g'(xk) 


9'k 


e'k 


9k 




-2.4 


-0.9399 


-0.4742 


-0.4955 


-0.7919 


0.1575 


-1.2 


1.0457 


5.5797 


4.3356 


0.8325 


0.2039 


0.0 


-0.7350 


5.4046 


-8.3534 


-0.5690 


0.2258 


1.2 


1.4159 


-2.6352 


2.8611 


0.5867 


0.5857 


2.4 


-1.0603 


-7.1265 


-5.7212 


-0.9327 


0.1203 


3.6 


0.2127 


-0.8092 


4.8046 


0.8036 


2.7786 



Table 7. Recovered samples of g' with two channels, r = 0.6, no noise 
added: exact samples (col. 2), computed samples and relative errors 
(col.s 3 and 4), computed samples and the relative errors (col.s 5 and 6) 
with regularization. 



Of course, by using a smaller oversampling parameter, one can reduce the condi- 
tion number of the matrix I—S , obtaining much better results: using r = 0.3 (hence 
t = 0.6), and the same set U = {-2,-1,0, 1,2,3}, one gets cond(I-S) = 1.92- 10 4 
and maximum absolute error on the samples of the signal equal to 0.0078. In Sec- 
tion 2 we have observed that, in the two-channel formula, the case r G (0, .5) implies 
oversampling in each channel separately, thus it can be used the one-channel for- 
mula separately for the signal and its derivative. However, the two alternative 
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are not numerically equivalent. To recover the same missing samples via the one- 
channel formula, we have taken r — 0.6 and U = {—2,-1,0,1,2,3} and we have 
obtained a maximum absolute error on the samples of the signal equal to 0.0274. 
This shows that the second alternative may be less efficient in case of consecutive 
samples. 




-25 -20 -15 -10 -5 5 10 15 20 25 



Figure 5. Original and reconstructed signal g, recover with regulariza- 
tion, for U = {-2,-1,0, 1,2,3} r = 0.6. 
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